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INTRODUCTION 

This is a final technical report of the project entitled "Crustal Deformation 
Along the San Andreas, California" supported by NASA Crustal Dynamics 
Program, for the period Feb. 1, 1991 - Jan. 31, 1992, based on funding from the 
Grant No. NAG 5-1398 awarded to the University of Michigan, Ann Arbor. 

The goal of this project is to achieve a better understanding of the regional 
and local deformation and crustal straining processes in western North America, 
particularly the effect of the San Andreas and nearby faults on the spatial and 
temporal crustal deformation behavior. Construction of theoretical models based 
on the mechanics of coupled elastic plate /viscoelastic foundation and large-scale 
crack mechanics provide a rational basis for the interpretation of seismic and 
aseismic anomalies and expedite efforts in forecasting the stability of plate 
boundary deformations. 

In the present period, special focus is placed on the three dimensional 
time-dependent surface deformation due to localized slippage in a . elastic layer 
coupled to a visco-elastic substrate. The numerical treatment of this problem is 
based on a 3-D boundary element technique derived in previous periods of 
research for elastic layered media. Extension to visco-elastic coupling demands 
the derivation of 3-D time-dependent Green's function. This has been 
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successfully achieved. This method has been applied to analyze the viscoelastic 
surface displacements due to a dislocated embedded patch. Surface uplift as a 
function of time and position are obtained. Comparisons between surface uplift 
for long and short dislocated patches are made. 


SUMMARY OF PROJECT ACHIEVEMENTS TODATE 


Major achievements for the duration of the project are summarized below: 


1 . Developed and completed a model of coupling between the elastic lithosphere 
and viscoelastic asthenosphere which incorporate the sub-mantle steady state 
motion as a long term driving force, and the transient loading and reloading of the 
plate boundary in earthquake cycles. This model suggests significant non- 
linearity in the stress-accumulation process over an earthquake cycle at a plate 
boundary. 

2. Developed and completed the modeling of crustal deformation (referred to as 
the Li-Rice model hereon) which varies with time and space at a strike-slip plate 
boundary embedded in the plate structure described in (1) above. 

3. From (1 ) and (2), and whole cycle San Andreas composite geodetic data, the 
surface elastic plate thickness is constrained to 20-30 km, and the viscoelastic 
relaxation time is constrained to 10-16 yr. 

4. The predicted surface velocity profiles are in good agreement with geodetic 
measurements at several locations along the San Andreas where such data is 
available. Locations studied include Point Reyes area and the Palmdale area. 
Predicted surface velocity profile at Palmdale area agree with the recently 
obtained field data. Note that model prediction is made independently of the 
geodetic field data indicated. 

5. Development of a model similar to the Li-Rice model for complex plate 
boundaries, including plate boundaries exhibiting surface fault creep, and plate 
boundaries with sub-parallel faults. Application of this model to the Parkfield area 
in Central California and to the Coachella-Valley area in Southern California is 
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performed. The predicted surface velocity for the Coachella-Valley area can be 
fitted ;to recent VLB I and ground based geodetic data. 

6. Comparison of the single line dislocation model of Savage and Burford with 
the Li-Rice model reveals that the locked depth and deep aseismic slip rate are 
significantly overestimated for the single line dislocation model using geodetic 
data as constraining information. Further, over the span of an earthquake cycle, 
the single fine dislocation model is shown to predict a shallower locked depth and 
a larger deep slip rate for a velocity profile associated with a time earlier in an 
earthquake cycle in comparison to that associated with a time later in an 
earthquake cycle. This is physically implausible. The suggestion from this 
analysis is that while the Savage/Burford single line-dislocation model is simple 
to use and understand, it also imposes severe limitation in direct physical 
interpretation of fault locked depth and deep aseismic slip rate, as well as other 
tectonic implications. 

Since the last technical report, our work has focused on the study of three 
dimensional effect of localized slippage in a elastic layer coupled to a visco- 
elastic substrate. This method has been applied to analyze the viscoelastic 
surface displacements due to a dislocated embedded patch. Surface uplift as a 
function of time and position are obtained. Comparisons between surface uplift 
for long and short dislocated patches are made. 

INVESTIGATION ON 3-D TIME DEPENDENT VISCOELASTIC RESPONSE 

The viscoelastic Green's functions considered were the displacement and stress 
fields due to a nucleus of strain located in an elastic plate (0<z<H) with a free 
surface at z=0 and a perfectly bonded lower (Figure. 1) viscoelastic halfspace at 
z>H. These Greens' functions will be used in surface integrals over slip patches 
of any arbitrary shapes, and can be used to deduce crustal deformations over a 
range of -time after the patch slippage. 

The viscoelastic model used for the lower halfspace is that of a Maxwell 
fluid described by 
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where Sjj and e,j are the deviatoric part of the stress and strain components Gjj 
and Ejj. A dot over a variable indicates rate with respect to time, p. and K are the 
shear and bulk modulus, and q s and rid are the shearing and volumetric 
viscosities respectively. For a quasi-static problem where the inertial terms are 
negligible, the only time dependence occurs in the boundary conditions and 
constitutive relations. It can be shown that the constitutive relations is similar in 
form in the Laplace transformed space to the elastic stress-strain relations in real 
space, but with the shear and bulk moduli replaced by moduli dependent on the 
Laplace variable "s". Thus the correspondence principle can be applied to solve 
the viscoelastic problem (Christensen, 1982). 


Within the context of the image method developed by us (Fares and Li, 
1988), and for the case when the nucleus of strain lies in the plate, field variables 
are only required in the plate, and where only the lower halfspace is viscoelastic, 
the only instances where moduli in Laplace transformed space which are 
dependent on the Laplace variable are in the "reflection" operators. It can be 
shown that there are in general four characteristic time scales for the response of 
a point source located in a plate perfectly bonded to one Maxwellian halfspace. 
The Laplace transform of a suddenly applied and then maintained unit excitation 
is 1/s. These considerations lead to insights into the behavior of the viscoelastic 
solutions. Specifically, the Laplace inverse of the terms in the solution after the 
correspondence principle is applied : 
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where -the lower case letter variables and Q are defined in terms of the moduli. 
Other terms are defined in Fares and Li (1988). After some algebra, the terms 
that are not multiplying the exponential functions in (3) - (5) can be simplified to 
the following 
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Therefore, the short time response of (3) - (5) (i.e. t -> 0+) corresponds to 
the halfspace being purely elastic and the long time response (from (6)-(8)) 
corresponds to a stress free boundary condition on the lower interface. 

The following observations are noted: i) The effective time scale of t m e' xt is 
longer than the effective time scale of t n e* xi when m>n, ii) the components of the 
product of the reflection matrices for "farther off" images contain polynomials in 
"t" with higher degrees multiplying the exponential terms, and iii) the relative 
contribution of "farther off" images increases as we get farther from the source 
location in the plate for fixed values of moduli. Thus the effective relaxation rate 
decreases as the distance from the source region increases. This concept could 
be associated with the stress diffusion phenomenon. 

As a numerical example of the above remarks, figures 2 and 3 
show the surface displacement components u x and u z respectively for embedded 
nuclei of strain (or moment sources) located at x=.y=0 and z=0.5H and having the 
normal vector n = (0,1,0) and (1 ,0,0) respectively and Displacement Discontinuity 
(DD) components Au = (1/pA). (1,0,0) and (1/pA). (0,0,1) respectively, where fi is 
the shear modulus and “A” is an “effective area”. The surface displacements are 
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plotted for y=z=0 and 0</H<4.0; the surface displacements are antisymmetric 
with respect to reflection of the x-axis. The time T at which the curves are 
plotted are cot = 0.0, 0:5, 1.0, 2.0, 3.0, 4.0, 50, 10.0, 15.0 and 30.0 for figure 2 
and, cot = 0.0, 1.0, 3.0, 5.0, 7.0, 10.0, 12.5 and 30.0 for figure 3; for both figures, 
the level of displacement increases with time (at least) in the range of 0<x/N<1 .0. 


The “diffusion" aspect of the displacement fields can be recognized with 
the help Of the following observations. In figure 2, consider two separate 
locations, -one relatively “close” to the source (e.g. at x/H=0.5) and one relatively 
“far” from the source (e g. at x/H=4.0). While the displacement seems to 
increase most at “early” times afterwhich this increase “slows down” for the closer 
location, the displacement increases relatively “faster” at “later” times for the 
“farther” location. The above remarks can be interpreted as being due to a 
“disturbance” which reaches and “saturates” the displacement field closer to the 
source and “diffuses" to reach farther off points from the source at “later" times. 
The time increments at which curves are plotted in figure 2 are not equal. 
Specifically, the last time step cot = 30.0 is chosen so that no significant change in 
the displacement magnitudes (for all locations 0<x/H<4.0) is noticed beyond that 
time step. The diffusion aspect in figure 3 is more readily apparent. A 
“disturbance” in the displacement magnitudes “seems” to be “penetrating" (or 
diffusing) farther from the source with increase in time. Again (for clarity of the 
figures), the time increments at which curves are plotted in figure 3 are not equal 
and cot = 30.0 is chosen so that no significant change in the displacement 
magnitudes (for all location 0<x/H<4.0) is noticed beyond that time step. 

INVESTIGATION OF VISCOELASTIC SURFACE DISPLACEMENTS DUE TO A 
DISLOCATED EMBEDDED PATCH 

A time-dependent model consisting of slippage (or dislocation) along a 
surface in an elastic. plate overlying a viscoelastic foundation is usually used in 
geophysical applications to study (or to take into account) postseismic surface 
deformations (i.e. surface deformations that occur "after" an earthquake). Two- 
dimensional "thrust" slippage and antiplane "strike slip" (both assume infinitely 
long faults) viscoelastic models have been repeatedly considered by many 
researchers (e.g. Nur and 'Mavko 1974, Smith 1974, Savage and Prescott 1978, 
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Thatcher and Rundle 1979, Lehner et al. 1981,'Li and Rice 1987 and Reilinger 
1986). 


In our research, the surface uplift due to uniform "thrust" slippage that 
occurs and is maintained along an inclined finite rectangular patch (two different 
lengths are considered:) in a plate has been analyzed. The geometry of the side 
and top views of those two patches are shown in figure 4. The aim of this study 
is to present surface .deformations due to a 3-D viscoelastic model of a layered 
earth, and to try and discern ;(if any.) deformation feature differences between a 2- 
D and a 3-D model. The short .patch represents a 3-D model and the long patch 
represents an approximate 2-D model. The calculation of the time-dependent 
surface uplift is obtained by the integration of the viscoelastic Green's functions 
discussed earlier over the surface of the patches using specialized integration 
techniques. 

The geometric parameters of the model are as shown in figure 4 (the short 
and long patches having lengths of 1 .2H and 5H respectively). Unit slippage (i.e. 
1H) is prescribed to occur uniformly over the surface of each patch; the surface 
uplift plots to be presented in this section (being linearly proportional to the 
amount of slippage along the patch) can be alternately viewed as being 
normalized by the slippage value. The elastic properties are such that the plate 
and the halfspace have identical shear modulus "G" and Poisson’ ratio "v"; the 
shear modulus "G" is taken to be "1" (arbitrary units) and Poisson's ratio "v" is 
taken to be 0.25. The units of the shear modulus does not affect the surface 
deformation values but serves to specify how much "force" has to be applied at 
the opposing slipping surfaces in order to produce the slippage (in geophysical 
applications, the shear modulus is said to affect the "moment" and "stress drop" 
on the patches). The volumetric deformation is purely elastic and the deviatoric 
deformation is a Maxwell model with a viscosity parameter rj = 20 (arbitrary 
units); the ratio of "r|/(2G)'' determines the characteristic time "t r " (time units 
depend on the units of both "n" and "G"). The above values of "ti/G” and "v" 
imply at = 2G/fi =0.1, Q = 0.5co, T2a = 0.33w, Qbi = 0.9co and Qb2 = 0.24co; the 
preceding values of Q, Qa, Qbi and ^B 2 show that the tour distinct exponential 
time-scales lie within an interval of "1/co" to "4.2/(o" for a model having the elastic 
properties of the plate and halfspace equal, and with Poisson’s ratio equal to 
0.25. 
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The surface uplift values along the X-axis (z=0, y=0 in figure 4) for 
-4.2<x/H<4.2 and up till sufficiently long limes (so that "steady state" has 
occurred in the range -0.7<x/H<0.7) are presented in figures 5 to 13 for both the 
short (figures 5 to 9) and long patch (figure -10) as well as for comparisons 
between the two patches (figures 1 1 to 13). A discussion of -these results will be 
next presented. 

First consider the uplift values due to the short patch. For convenience 
the deformation pattern will be divided into three sections. These sections 
numbered 1 to 3 lie roughly between -4.2<x/H<0.6, -0.6<x/H<o.04 and 
0.4<x/H<4.2 respectively. The material lying at x/H<0.3 is being "pushed” onto a 
narrow "gap" between the free surface and the closest edge of the rectangular 
patch to that free surface. This "squeezing" effect causes a sharp peak in the 
uplift values to occur in section "2". However, the material lying at x/H<0.3 is 
being "anchored" to a halfspace which is steadily relaxing, and hence with time, 
the uplift values in section "2" move in the direction of "imposed slippage" on the 
upper surface of the patch. Finally, material is being driven below (or being 
pushed under) the patch (Starting from x<0.3 to x>-0.3) and this material is being 
resisted with a "floor" which is steadily relaxing; hence the uplift values in section 
"3" move in the direction of "imposed slippage" on the lower surface of the patch. 
The evolution of uplift with time is quite complex in detail, as can be seen in 
figures 6-9. Some features that can be discerned in those figures is an 
exponential decay of surface deformation towards a "steady state" in section "1" 
(figures 7 and 8), as well as some diffusional behavior which can be deduced 
from the "more elongated" and "longer lived activity" time behavior of locations 
that are farther off from the patch (figures 6 and 9); however, from the 
convergence studies of chapter 2, the numerical accuracy of the surface 
deformations of "faraway" locations are not very certain, and exact judgment on 
the results cannot be made. The next discussion will compare the "short" (3-D) 
versus the "long" (2-D) surface uplifts. 

The surface uplift for the "long" patch is shown in figure 10, and is not 
qualitatively much different from the surface uplift due to the "short" patch. A 
general observation is that the surface uplift values for the "long” patch has 
higher peaks and varies more with time than the "short" patch; this observation is 
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probably due to the higher compliance inherent in a longer patch (by having less 
"side supports"). Figures 1 1 to 1 3 show comparisons at short, intermediate and 
long time values between the uplift values of the "short" and "long" patch. Tihe 
higher "inherent compliance" of the longer patch causes the "anchoring" in 
section "1" (of the plots) to be less effective and the "driving under" in section "3" 
(of the plots) to be more effective than for the shorter patch. 
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Surface Uplift vs Position (Short Patch) 
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Figure 5 
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Surface Uplift vs Time 
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Surface Uplift vs Time 
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Figure 7 
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Surface Uplift vs Position (Long Patch) 
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Figure 10 
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Surface Uplift (t/tr = 0.0) 
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Figure 1 1 
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Surface Uplift (t/tr = 5.0) 



Figure 12 
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Surface Uplift ( t/tr = 30.0) 
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Figure 13 
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